{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "n = 60  m = 400\n"
     ]
    }
   ],
   "source": [
    "# import necessary modules\n",
    "import numpy as np               \n",
    "import matplotlib.pyplot as plt\n",
    "from IPython.display import display\n",
    "import os\n",
    "\n",
    "pwd = os.getcwd()\n",
    "\n",
    "# write your OWN PC folder path for fdir\n",
    "FOLDER = 'VESSEL_SED_CHAN' # WRITE NAME OF PROJECT FOLDER\n",
    "fdir = os.path.join(pwd,FOLDER,'output')\n",
    "\n",
    "# Getting depth file and determining domain dimensions\n",
    "\n",
    "depFile = os.path.join(fdir,'dep.out')\n",
    "dep = np.loadtxt(depFile)\n",
    "[n,m] = dep.shape\n",
    "print('n =',n, ' m =',m)\n",
    "\n",
    "N = 2*n+1\n",
    "M = m\n",
    "\n",
    "# discretization\n",
    "dx = 1.0\n",
    "dy = 1.0\n",
    "\n",
    "# x and y field vectors \n",
    "x = np.arange(0,M*dx,dx)\n",
    "y = np.arange(0,N*dy,dy)\n",
    "\n",
    "files=[20]\n",
    "time=['40']\n",
    "\n",
    "ETA = np.zeros([N, M])\n",
    "CH = np.zeros([N, M])\n",
    "BEDS = np.zeros([N, M])\n",
    "BEDB = np.zeros([N, M])\n",
    "BB = np.zeros([N, M])\n",
    "DEP = np.zeros([N, M])\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/rdchlgs8/anaconda/lib/python3.6/site-packages/ipykernel_launcher.py:46: RuntimeWarning: divide by zero encountered in log10\n",
      "/Users/rdchlgs8/anaconda/lib/python3.6/site-packages/matplotlib/colors.py:496: RuntimeWarning: invalid value encountered in less\n",
      "  cbook._putmask(xa, xa < 0.0, -1)\n"
     ]
    }
   ],
   "source": [
    "## Plot Output\n",
    "size = 15\n",
    "\n",
    "\n",
    "for num in range(len(files)):\n",
    "    fnum = '%.5d' % files[num]\n",
    "\n",
    "    eta = np.loadtxt(os.path.join(fdir,'eta_'+fnum))\n",
    "    mask = np.loadtxt(os.path.join(fdir,'mask_'+fnum))\n",
    "    ch = np.loadtxt(os.path.join(fdir,'C_'+fnum))\n",
    "    beds = np.loadtxt(os.path.join(fdir,'DchgS_'+fnum))\n",
    "    bedb = np.loadtxt(os.path.join(fdir,'DchgB_'+fnum))\n",
    "    \n",
    "    # Removing masked regions from plot\n",
    "    eta[np.where( mask < 1)] = np.nan\n",
    "    ch[np.where( mask < 1)] = np.nan\n",
    "    \n",
    "    ETA[0:n,:] = eta[:,:]\n",
    "    ETA[n:-1,:] = eta[::-1,:]\n",
    "    CH[0:n,:] = ch[:,:]\n",
    "    CH[n:-1,:] = ch[::-1,:]\n",
    "    \n",
    "    BEDS[0:n,:] = beds[:,:]\n",
    "    BEDS[n:-1,:] = beds[::-1,:]\n",
    "    BEDB[0:n,:] = bedb[:,:]\n",
    "    BEDB[n:-1,:] = bedb[::-1,:]\n",
    "    \n",
    "    BB = BEDS+BEDB\n",
    "    \n",
    "    # plot eta (surface elevation)\n",
    "    fig = plt.figure(figsize=(16, 12))\n",
    "    \n",
    "    ax1 = fig.add_subplot(3,1,1)\n",
    "    ax1.tick_params(axis='both', which='major', labelsize=size)\n",
    "    ax1.tick_params(axis='both', which='minor', labelsize=size)\n",
    "    plt.pcolor(x,y,ETA, cmap='jet')\n",
    "    ax1.axis([0, 400, 0, 120])\n",
    "    h_bar = plt.colorbar(label=r'$\\eta (m)$')\n",
    "    plt.clim(-2.0, 2.0)\n",
    "    plt.ylabel('y (m)',fontsize=size)\n",
    "    \n",
    "    # plot log10 \n",
    "    ax2 = fig.add_subplot(3,1,2)\n",
    "    ax2.tick_params(axis='both', which='major', labelsize=size)\n",
    "    ax2.tick_params(axis='both', which='minor', labelsize=size)\n",
    "    plt.pcolor(x,y,np.log10(CH), cmap='jet')\n",
    "    ax2.axis([0, 400, 0, 120])\n",
    "    l_bar = plt.colorbar(label='log10(C) (mg/L)')\n",
    "    plt.clim(0, 0.5)\n",
    "    plt.ylabel('y (m)',fontsize=size)\n",
    "    \n",
    "    # plot BB \n",
    "    ax3 = fig.add_subplot(3,1,3)\n",
    "    ax3.tick_params(axis='both', which='major', labelsize=size)\n",
    "    ax3.tick_params(axis='both', which='minor', labelsize=size)\n",
    "    plt.contourf(x,y,BB,10, cmap='jet')\n",
    "    ax3.axis([0, 400, 0, 120])\n",
    "    bb_bar = plt.colorbar(label='dZ (m)')\n",
    "    plt.clim(-0.002, 0.002)\n",
    "    plt.ylabel('y (m)',fontsize=size)\n",
    "    plt.xlabel('x (m)',fontsize=size)\n",
    "    \n",
    "    plt.savefig(os.path.join(pwd,'vessel_eta_c_dz_%s.jpg'%(time[num])), dpi=400) \n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "[nn,mm]=np.shape(BB)\n",
    "\n",
    "DEP[0:n,:]=dep[:,:]\n",
    "DEP[n:-1,:]=dep[::-1,:]\n",
    "\n",
    "B = np.zeros(np.shape(BB[0,:]))\n",
    "           \n",
    "# plot average change in morpho\n",
    "for j in range(nn):\n",
    "    B[j] = np.mean(BB[j,:])\n",
    "\n",
    "fig = plt.figure(figsize=(16, 10))\n",
    "ax1 = fig.add_subplot(2,1,1)\n",
    "plt.plot(x,B,linewidth=2.0)\n",
    "ax1.grid(linestyle='--', linewidth=1) \n",
    "plt.xlabel('y (m)')\n",
    "plt.ylabel('Averaged Bed Change (m)')\n",
    "ax1.axis([0 ,120, -0.004, 0.002])\n",
    "\n",
    "ax2 = fig.add_subplot(2,1,2)\n",
    "plt.plot(y,-1*DEP[:,0],linewidth=2.0)\n",
    "plt.plot([0, 120],[0, 0],'b--')\n",
    "ax2.grid(linestyle='--', linewidth=1) \n",
    "plt.xlabel('y (m)')\n",
    "plt.ylabel('Initial Depth(m)')\n",
    "ax2.axis([0, 120, -5.2, 2])\n",
    "\n",
    "plt.savefig(os.path.join(pwd,'vessel_bedchg_%s.jpg'%(time[num])), dpi=400)\n",
    "             \n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "hide_input": false,
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.1"
  },
  "nbTranslate": {
   "displayLangs": [
    "*"
   ],
   "hotkey": "alt-t",
   "langInMainMenu": true,
   "sourceLang": "en",
   "targetLang": "fr",
   "useGoogleTranslate": true
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
